Skip to content

Issue 525 FE Basis refactor implementation#561

Open
zasexton wants to merge 116 commits into
SimVascular:mainfrom
zasexton:issue-525
Open

Issue 525 FE Basis refactor implementation#561
zasexton wants to merge 116 commits into
SimVascular:mainfrom
zasexton:issue-525

Conversation

@zasexton

@zasexton zasexton commented Jun 2, 2026

Copy link
Copy Markdown
Collaborator

Current situation

Tracks #525.

The solver currently relies on legacy table-driven shape-function paths in nn.cpp, which makes basis evaluation, Hessian support, node-ordering validation, and parity testing difficult to extend. This PR introduces a self-contained FE Basis layer while preserving the existing solver-facing storage contracts.

Release Notes

  • Added a new FE Basis module with Lagrange and Serendipity basis support for mapped solver elements, including values, gradients, Hessians, node-ordering conventions, factory creation, and cache support.
  • Added FE math and quadrature utilities used by the new basis implementation.
  • Routed supported nn::get_gnn and nn::get_gn_nxx volume/face evaluations through the new FE Basis adapter.
  • Added SVMP_BASIS_MODE=auto|legacy|fe to compare legacy and FE Basis evaluation paths at runtime.
  • Added typed FE/basis exception paths for unsupported elements, invalid configurations, and backend failures.
  • Added a basis comparison harness for running pytest integration cases under legacy and FE modes and comparing outputs, timing, fields, and iteration data.
  • Build migration: solver code now requires a C++20-capable compiler.

Documentation

  • Added Doxygen-style documentation to the new FE Basis, FE Math, and Quadrature APIs.
  • No solver XML migration is expected for users; basis mode selection is controlled through SVMP_BASIS_MODE.

Testing

  • Added GoogleTest coverage under tests/unitTests/FE/Basis for Lagrange/Serendipity basis evaluation, Hessians, cache/factory behavior, error paths, node ordering, solver adapter parity, and supported mapped element coverage.
  • Added GoogleTest coverage under tests/unitTests/FE/Math for matrix/vector operations, expression helpers, math constants, and dense linear algebra.

Code of Conduct & Contributing Guidelines

@zasexton zasexton self-assigned this Jun 2, 2026
@zasexton

zasexton commented Jun 3, 2026

Copy link
Copy Markdown
Collaborator Author

From the CI test cases we got the following:

  • macOS: 12 failed, 227 passed
  • Ubuntu: 21 failed, 218 passed
  • Failures are all FSI / FSI-ustruct pipe cases.

The affected meshes are HEX8. Current code allocated Nxx but did not populate HEX8 second derivatives; this branch now computes nonzero HEX8 Hessians via FE Basis. FSI uses those in Code/Source/solver/fsi.cpp:191, so the stored legacy reference VTUs no longer match.

@zasexton

zasexton commented Jun 3, 2026

Copy link
Copy Markdown
Collaborator Author

In the FSI fluid assembly, Nxx is passed into the element kernel and used in terms that contribute to the assembled equations, especially VMS/stabilization terms.

  • Code/Source/solver/fluid.cpp:1841: Nwxx is used to build uxx at lines 1841-
    1860.

  • Code/Source/solver/fluid.cpp:1880: uxx feeds d2u2.

  • Code/Source/solver/fluid.cpp:1933: uxx also feeds es_x, then mu_x at lines
    1948-1963.

  • Code/Source/solver/fluid.cpp:2027: d2u2 and mu_x enter rS.

  • Code/Source/solver/fluid.cpp:2037: rS enters up.

  • Code/Source/solver/fluid.cpp:2103: up-dependent quantities are used when
    assembling lR, the local residual, at lines 2103-2106.

That means:

  • in the current code, it solves a discrete problem where HEX8 second-derivative contributions are effectively zero or absent when using the old basis functions

  • when using the new FE Basis functions, it solves a different discrete problem where HEX8 mapped Hessian contributions are present

Therefore the Hessian is not just participating in the tangent assembly and may not necessarily converge on the same solution as the current reference vtu. This is also why we do not see such CI failures for test cases that use Tet elements because TET4 have zero reference second derivatives or for TET10 Hessian values are provided. For such elements, introduction of the new Basis function infrastructure did not cause a change in the residual being solved for these cases.

@zasexton

zasexton commented Jun 3, 2026

Copy link
Copy Markdown
Collaborator Author

Overall, we see that providing the non-zero Hessian values slightly increases the number of linear iterations required but does not impact the number of nonlinear iterations for the FSI test cases considered. More importantly, this difference in how the Hessian is used in the formulated residual suggested that the problem being solved is slightly different depending on if a user is providing brick elements versus TET10 (which shouldn't be the case)

convergence-iteration-behavior

@zasexton

zasexton commented Jun 4, 2026

Copy link
Copy Markdown
Collaborator Author

Confirmed. The difference in accuracy for the FSI test cases when employing the new FE Basis functions is due to the non-zero Hessian values for the reference elements. We can see that when we zero-out the values for the Hessian on the HEX8 cases that we recover the same behavior as the current solver. @ktbolt do you think it would make sense to update the reference solutions for these test cases or is there an alternative strategy that we should look at?

hex8-hessian-confirmation-by-case

@aabrown100-git

Copy link
Copy Markdown
Collaborator

@zasexton Were Hessians also zeroed out in the Fortran svFSI?

I think we should update reference solutions unless there is a good reason to use zero Hessians. @sujaldave, @hanzhao2020 Is there any reason (in VMS or something) to force $u_{xx}$ to be zero in HEX8 or higher order elements?

@zasexton

zasexton commented Jun 4, 2026

Copy link
Copy Markdown
Collaborator Author

@aabrown100-git yes, the Hessian was also effectively zero for the HEX8 in the original fortran svFSI. Only TET10, TRI6, QUD8, QUD9, and LIN2 had cases in the GETGNNxx subroutine. All other cases return the zero-initialized Hessian.

Regenerate affected FSI and FSI-ustruct HEX8 result_005.vtu references for the FE Basis path with nonzero HEX8 Hessian contributions.

Update the pipe_3d PETSc and Trilinos references to match the base pipe_3d reference, preserving the existing shared-reference pattern across linear algebra variants.
@codecov

codecov Bot commented Jun 8, 2026

Copy link
Copy Markdown

Codecov Report

❌ Patch coverage is 92.01371% with 303 lines in your changes missing coverage. Please review.
✅ Project coverage is 72.98%. Comparing base (9632608) to head (49ba342).

Files with missing lines Patch % Lines
Code/Source/solver/nn.cpp 66.17% 69 Missing ⚠️
Code/Source/solver/Parameters.cpp 38.94% 58 Missing ⚠️
tests/unitTests/FE/Basis/test_BasisErrorPaths.cpp 81.25% 48 Missing ⚠️
tests/unitTests/FE/Basis/test_LagrangeBasis.cpp 94.19% 28 Missing ⚠️
tests/unitTests/FE/Basis/test_SerendipityBasis.cpp 97.24% 19 Missing ⚠️
tests/unitTests/FE/Basis/test_BasisHessians.cpp 93.33% 15 Missing ⚠️
Code/Source/solver/fs.cpp 37.50% 10 Missing ⚠️
...Source/solver/FE/Basis/NodeOrderingConventions.cpp 98.52% 7 Missing ⚠️
Code/Source/solver/FE/Math/DenseLinearAlgebra.cpp 96.29% 7 Missing ⚠️
Code/Source/solver/FE/Basis/LagrangeBasis.cpp 98.21% 6 Missing ⚠️
... and 12 more
Additional details and impacted files
@@            Coverage Diff             @@
##             main     #561      +/-   ##
==========================================
+ Coverage   69.44%   72.98%   +3.54%     
==========================================
  Files         181      237      +56     
  Lines       34072    37670    +3598     
  Branches     5930     6539     +609     
==========================================
+ Hits        23662    27495    +3833     
+ Misses      10273     9940     -333     
- Partials      137      235      +98     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

@michelebucelli michelebucelli left a comment

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you @zasexton, this is lots of work! I added some other comments.

Also good job with the documentation, it's very detailed.

Comment thread Code/Source/solver/FE/Basis/BasisFactory.h Outdated
Comment thread Code/Source/solver/FE/Basis/BasisFactory.h Outdated
Comment thread Code/Source/solver/Core/Exception.h
Comment thread Code/Source/solver/FE/Basis/BasisExceptions.h
Comment thread Code/Source/solver/FE/Basis/BasisExceptions.h Outdated
Comment thread Code/Source/solver/FE/Basis/LagrangeBasis.cpp Outdated
Comment thread Code/Source/solver/FE/Common/Types.h
Comment thread Code/Source/solver/FE/Common/Types.h Outdated
Comment thread Code/Source/solver/FE/Common/Types.h Outdated
Comment thread Code/Source/solver/nn.cpp Outdated
zasexton added 21 commits June 25, 2026 14:58
- include the basis family in the documented basis identity (topology, order, family)
- replace the order() normalization wording with a plain description and example
- add cross-reference links to the module entities in the group documentation
- add topology and family examples to the basis object contract
- refer to shape function values in the BasisFunction class description
Add an explanatory sentence and a concrete example to each derived basis
exception, drawn from its real throw sites, so the specific types guide
developers on when to raise or expect them.
…ation history

Name the orders and elements that use equispaced nodes directly, instead of
referring to what previous element layouts kept, so the note stays accurate as
the implementation evolves.
create() and create_default_for() now return std::unique_ptr<BasisFunction>,
which has simpler ownership semantics and is cheaper than shared_ptr. The single
caching consumer keeps a shared_ptr cache; the returned unique_ptr converts to it
on insertion. All FE/Basis unit tests pass.
…erloads in the base

The span *_to evaluators are now the primitives a concrete basis implements:
evaluate_values_to is pure virtual, evaluate_gradients_to/evaluate_hessians_to
default to reporting not-implemented, and a protected evaluate_all_to provides the
single-pass combined evaluation. The std::vector overloads (evaluate_values,
evaluate_gradients, evaluate_hessians, evaluate_all) are implemented once on the
base class, so LagrangeBasis and SerendipityBasis no longer duplicate those
wrappers.

This inverts the previous arrangement where the vector form was the primitive and
the base supplied an allocate-and-copy span fallback. A minimal basis now
overrides the span primitive instead of the vector form. The BasisFunction test
helpers are updated accordingly, and the former fallback test is rewritten to
verify the base vector overloads forward to the span primitives. All FE/Basis
unit tests pass.
…tation

- note that the vector evaluators use an output argument so callers can reuse a
  buffer across evaluations rather than allocating per call
- document that the BasisTraits classifiers return -1 / Unknown as sentinels that
  callers validate into exceptions (the constexpr noexcept helpers cannot throw)
- keep the BasisRequest::topology field-ordering note out of the rendered docs
…roup

Replace the @cond INTERNAL exclusion with a documented internal group whose
@warning states the declarations are internal (use basis_factory and
BasisFunction::nodes() instead) and may change. Core developers now get the
rendered documentation while model-level callers are clearly steered away.
… Vector

Document that FE/Math/Vector.h is a fixed-size, compile-time-length element-level
vector in svmp::FE::math, distinct from and not a replacement for the legacy
dynamic global ::Vector in solver/Vector.h.
…nd ElementType bands

- explain that the Mesh library is an optional external module and FE builds
  standalone with fallback types when it is absent
- record why GlobalIndex stays a plain alias (raw PETSc/Trilinos interop) and that
  DofIndex is the strong wrapper
- document that the explicit ElementType values are intentional bands with a
  uint8_t Unknown sentinel
…e_label

Add Doxygen to the dense linear-algebra functions and the DenseLUSolver type, and
rename the label argument (and the DenseLUSolver member) to error_message_label so
its role as an error-message prefix is self-documenting. All FE math unit tests pass.
Add a default case that raises BasisElementCompatibilityException for a solver
element type with no FE mapping, so a missing mapping fails loudly instead of
relying on the unhandled-enum compiler warning. The deliberate NA/PNT/NRB cases
still return std::nullopt.
MathJax 2 was pinned as a default when the FE documentation was added and is now
end-of-life. Move to the maintained MathJax 3 (chtml output, mathjax@3 CDN, ams
extension); doxygen emits the v3 bootstrap and the AMS math used in the FE docs is
supported.
…mbers

Generating HTML (not just parsing) surfaced reference/doc issues the no-output
checks missed:
- qualify the @ref class/struct/enum links in the Basis group doc so they resolve
- use plain text for the basis_factory namespace mention (no resolvable @ref)
- avoid the ::Vector explicit-link request in the Math vector doc
- add @return tags and document the NodeOrderingConventions entities that became
  visible when the @cond INTERNAL exclusion was removed
Net FE Doxygen warnings drop from 47 to 19 (the rest pre-existing).
…assertion helpers

- add svmp::NotImplementedException and svmp::IndexOutOfRangeException (CoreException-derived)
- default the ExceptionT template parameter of not_implemented() and check_index()
  to those types, so they can be used without naming an exception while existing
  explicit callers are unaffected
- document the helpers, including that check vs check_arg are the same check named
  for intent and that throw_if is the logical inverse of check, and when to use
  each not_implemented overload
Note: svmp::FE::NotImplementedException stays for FE code needing FEException
ancestry. All FE/exception-helper unit tests pass.
@zasexton zasexton marked this pull request as ready for review June 30, 2026 01:53

@claude claude Bot left a comment

Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Claude Code Review

This pull request is from a fork — automated review is disabled. A repository maintainer can comment @claude review to run a one-time review.

@zasexton

Copy link
Copy Markdown
Collaborator Author

I believe we should be ready to merge this PR into the repo. I believe the basis refactor is complete and the documentation is sufficient for future users and developers. the code has good coverage via unit testing and passes all current integration tests.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

OOP Refactor Object-Oriented Programming Refactor of Code

Projects

None yet

Development

Successfully merging this pull request may close these issues.

Implementing Basis Module

4 participants